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I.  Introduction 

The  Thaerocomp  Technical  Corporation  (TTC)  has  proposed  to  develop  next-generation  software 
that  will  revolutionize  practical  CFD-based  analyses  by  providing  the  capability  for  seamless 
multidisciplinary  simulation  using  highly  accurate  procedures.  The  software  identifies, 
integrates,  and  then  systematically  upgrades  common  elements  of  major  research  codes,  to  yield  a 
versatile  high-order  (Lele,  1991)  software  suite  that  can  be  used  to  simulate,  in  a  coupled  fashion, 
various  elements  of  turbulence,  electromagnetics,  aeroacoustics,  magnetogasdynamics,  fluid- 
structure  interaction,  and  chemically  reacting  flows.  An  elegant  and  intuitive  graphical  user 
interface  (GUI)  will  be  developed  by  leveraging  experience  derived  from  our  flagship  product 
INSTED.  TTC’s  new,  tightly  integrated  and  versatile  product  will  reduce  the  simulation  time  by 
orders  of  magnitude  and  drastically  reduce  the  design  cost  of  modern  engineering  systems 
requiring  multidisciplinary  optimization.  The  proposed  software  will  facilitate  out-of-the-box 
approaches  by  providing  a  validated  realistic  virtual  environment,  thus  impacting  all  major  thrusts 
of  the  Air  Vehicles  Directorate,  including  access-to-space,  sustainment,  and  unmanned  air 
vehicles. 

The  present  document  summarizes  the  accomplishments  of  the  Phase  I  effort,  which  successfully 
and  firmly  established  the  feasibility  of  the  concept  proposed  by  TTC.  In  specifics,  the  following 
tasks  were  successfully  carried  out  in  Phase  I: 

•  Identification  of  several  existing,  but  not  necessarily  high  order,  upgradable  target  codes,  with 
a  synopsis  of  the  numerical  characteristics  of  the  codes.  For  each  code,  several  phone  and 
written  contacts  were  made  in  order  to  obtain  the  current  status  of  the  code. 

•  Selection  of  the  base  code  for  the  development  of  the  high  order  multidisciplinary  software. 
The  specific  criteria  for  selection  include  a)  the  type  of  numerical  method  (finite 
element/volume,  finite  difference,  spectral  methods),  b)  capabilities  of  the  codes,  c)  grid 
topology  (structured/unstructured),  d)  the  structure  of  the  code,  vis-a-vis  the  incorporation  of 
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chimera/overset  capabilities  for  competitive  modeling  of  complex  geometry  and  relative 
motion,  e)  the  availability  of  the  code. 

•  Identification  of  the  numerical  kernels  in  the  codes  and  development  of  a  generalized  flux 
function  that  supports  all  the  applications. 

•  Implementation  of  an  elegant  and  intuitive  graphics  user  interface  (GUI)  to  demonstrate  the 
feasibility  of  a  user-friendly  interface  for  the  high  order  multidisciplinary  software,  which  is 
an  absolute  necessity  for  commercialization  of  the  final  product. 

•  Development  and  testing  of  a  preliminary  version  of  the  proposed  software.  A  demonstration 
for  MGD  and  Aeroacoustics  is  presented  in  terms  of  the  computational  efficiency  and 
accuracy  of  computation  when  compared  with  codes  that  are  dedicated  to  these  applications. 

•  Investigation  of  marketing  efforts  for  the  commercial  sector  during  Phase  I,  which  shows  that 
considerable  interest  exists  for  the  proposed  product  because  of  its  potential  to  ensure  fast 
time-to-market. 

II:  Identification  of  Target  Codes 

The  target  codes  for  upgrade  have  been  identified,  as  were  the  kernels  in  these  codes.  Some  of 
these  codes  are  listed  in  the  table  below,  together  with  a  synopsis  of  their  numerical 
characteristics.  Several  contacts  were  made  to  find  out  the  current  state  of  the  various  codes.  Dr. 
Pieter  Buning  (NASA  Langley)  provided  an  overview,  current  status,  and  the  relationship 
between  ARC3D  code  (originally  developed  by  Tom  Pulliam),  XAIR  (derivative  of  ARC3D), 
OVERFLOW  (also  a  derivative  of  ARC3D),  CFL3D  (developed  from  NASA  Langley),  and 
WIND  (originally  developed  by  Douglas  Ray  Cosner/McDonnell  Douglas).  Greg  Power 
(Contractor  at  AEDC)  provided  detailed  information  to  TTC  on  the  current  status  of  WIND, 
which  is  the  main  production  flow  solver  for  the  NPARC  Alliance. 


ARC3D 

FDL3DI 

XAIR 

OVERFLOW 

WIND  (NPARC) 

(NASA) 

(WPAFB) 

(AEDC) 

(NASA) 

(Boeing) 

Method 

FD 

FD 

FD 

FD 

FD  (Node-based 

FV) 

Spatial 

Scheme 

Upwind, 

centered 

Compact, 

centered, 

upwind 

Upwind, 

centered 

Upwind,  centered 

Upwind,  centered, 

Spatial 

accuracy 

2”'',  probably 
(Roe) 

2”“,  4'^  6‘^ 
etc. 

2"‘*,  probably 
3^‘*  (Roe) 

probably  3*^^ 
(Roe) 

probably  3*^^ 
(Roe) 

Temporal 

scheme 

Beam- 

Warming 

Beam- 

Warming, 

RK4 

Partial 
factorization 
+  Gauss- 
Seidel 

Beam- Warming 

Partial  factorization 
+  point  Jacobi  or 
Gauss-Seidel 

Time 

Accuracy 

Steady  state, 
second  order 
with  sub¬ 
iterations 

RK4,  2"“ , 

3'^‘* ,  time  - 
accurate 

2''^ ,  time- 
accurate 
moving 
grids 

Steady  state, 

2“** -order  with 
sub-iterations 

First-order, 

essentially 

Stabilization 

Explicit  4*  /2"'‘ 
dissipation, 
Inherent  Roe 

Compact 
filters, 
explicit  4/2 

HLLE  (Roe¬ 
like) 

upwinding 

Explicit  4*  IT^ 
dissipation, 
inherent  Roe 

Upwinding 

Multiblock 

No 

Yes 

Yes 

Yes 

Yes 

Chimera, 

Overset 

Capabilities 

None 

Limited 

Extensive 

Extensive 

Extensive 
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Turbulence 

modeling 

RANS 

LES 

RANS 

RANS 

RANS 

Matrix 

method 

Scalarized 

pentadiagonal 

TDMA 

N/A 

Scalarized 

pentadiagonal 

N/A 

Mesh 

interpolation 

tools 

Limited 

Limited 

Extensive 

Extensive 

Extensive 

Extent  of 
dissipation 

Significant 

Slight 

Significant 

Significant 

Significant 

Applications 

Standard 
Navier-Stokes, 
research, 
mostly  in- 
house 

Acoustics, 
DNS/LES, 
CEM,  MHD, 
mostly  in- 
house 

In-house 
only,  not 
supported 

Extensive  industry 
use 

Part  of  NPARC, 
extensive  industry 
use 

On  the  basis  of  the  numerical  characteristics  of  the  reviewed  code,  vis-a-vis  the  goal  of  this 
project,  we  selected  the  FDL3D  code  (Gaitonde  and  Visbal,  1998)  as  the  basis  for  the 
development  of  the  multidisciplinary  CFD  tool.  The  specific  criteria  used  include  the  type  of 
numerical  method,  the  capabilities  of  the  target  code,  grid  topology  in  the  code,  chimera/overset 
capability  and  the  availability  of  the  code.  Note  that  none  of  the  versions  of  the  FDL3D  code  is 
capable  of  multidisciplinary  simulation.  They  also  do  not  currently  provide  high-order  shock¬ 
capturing  capabilities.  Furthermore,  the  overset/chimera  capability  is  not  as  advanced  as  that  in 
some  of  the  other  codes  that  we  reviewed,  notably  OVERFLOW  and  WIND.  The  types  of 
boundary  conditions  currently  supported  in  FDL3D  are  limited,  for  example,  when  compared 
with  WIND.  These  4  areas  were  automatically  added  to  the  technical  focus  for  the  SBIR  project. 
The  WENO  scheme  (Shu,  1997)  was  selected  for  high  order  shock  capturing. 

Ill:  Governing  Equations  and  mathematical  model 

The  Phase  I  proposal  promised  to  develop  a  preliminary  generalized  framework  for  the 
multidisciplinary  software,  including  a  generalized  flux  function,  and  to  select  the  platform  for 
the  high  order  calculations  and  for  stabilization  in  flow  fields  with  shocks  and/or  strong 
discontinuity.  It  was  also  planned  to  test  a  preliminary  version  of  the  software  suite. 

A  preliminary  version  of  the  proposed  module  that  supports  generalized  flux  functions  and  source 
terms  was  developed  in  Phase  I.  The  module  uses  the  same  high  order  procedure  for  evaluating 
the  derivatives  for  the  various  flux  functions  (source  terms).  The  analysis  of  the  efficiency  and 
accuracy  of  the  generalized  flux  procedure  relative  that  of  the  individual  single  applications  was 
also  started  in  Phase  I.  Each  of  the  supported  applications,  in  curvilinear  coordinates,  can  be 
written  in  the  form 


dx^  dF  dff 

dr  St] 


dC  d^  d7] 


+  ^  +  S' 


where  x  is  the  time  domain,  (^,  rj,  Q  are  the  transformed  coordinates  that  respectively  correspond 
to  the  physical  coordinates  (x,y,z).  F\G  \  and  H'  are  the  inviscid  fluxes,  while  F\,G\,  and  H’  ^ 
are  the  fluxes  associated  with  the  stresses  for  a  particular  application.  S'  is  3.  generic  source  term. 
Some  manipulations  are  sometimes  required  in  order  convert  the  stresses  into  the  foregoing  form, 
notable  examples  being  the  MaxwelVs  stresses  in  MGD.  The  various  flux  functions  that  need  to 
be  “combined”  are  listed  below,  using  F'  as  an  example.  The  generalized  flux  functions  for  some 
of  the  multidisciplinary  applications  are  presented  subsequently. 
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DNS,  non-MHD\  When  a  handle  to  the  module  indicates  direct  numerical  simulation,  or  the  fluid 
component  of  a  fluid-structure  interaction,  the  module  assumes  the  following  definitions 


p 

(  \ 

pu 

puU  + 

pv 

pvU  +  ^^p 

pw 

pmU+^jp 

pE 

\  p{E+p)u-i,p  y 

where 


LES,  non-MHD\  For  non-MHD,  LES  problems,  the  following  definitions  apply  (Germano  et  ah, 
1991;  Rizetta  et  al,  2000,  2001): 


P 

(  \ 

pu 

r  = 

pm 

pwU  + 

pE 

^  pEU  +  ^^Uip  j 

where 


MGD:  For  a  simplified  MGD  model,  i.  e.,  assuming  chemical  and  thermal  equilibrium,  the 
appropriate  definitions  are  (Gaitonde,  1999): 


/  pu 


9 

pu 

pvU  + 

pw 

pE 

E  =  - 
'  J 

pwU  + 

{pE+v)U-i^R^{\5  ■  f-)B^-i  R^{i3  .  f-)B^ 

4  -^^y)  +  L 

y 

i^iuBy  -vB^)  +  i^(wBy  -vB^) 

\  i^i'>jE^-wB^)+4y(vB^-vjBy) 

\ 


J 


The  standard  notations  apply  (see  reference). 

Species  conservation:  When  reacting  or  non-reacting  species  have  to  be  simulated,  the  module 
must  branch  to  the  following  definitions: 


’  P'Pl  ’ 

(  P^iU  \ 

1 

P^2 

P'p2^ 

^'  =  7 

_  P4>r.  _ 

,  r  = 

(  pKU  / 
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where,  as  before,  we  can  define 


^  +  + 

Note  that  the  species  fluxes  could  easily  be  appended  to  any  of  the  foregoing  definitions,  in  a 
coupled  approach,  as  opposed  to  the  segregated  manner  implied  here. 

CEM:  For  the  time-dependent  electromagnetics,  the  modules  sends  controls  to  the  Maxwell's 
equations  for  the  electromagnetic  field  in  free-space,  which  can  be  written  as  (Shang  and 
Gaitonde,  1995,  Gaitonde  and  Shang,  1997): 


0 

0 

1 

0 

= 

J 

J. 

J. 

In  these  expressions,  B  =  (Bx,  By,  Bz)  is  the  vector  of  magnetic  field  density,  D  =  (Dx,  Dy,  Dz)  is 
the  vector  of  electric  displacement,  pm  is  magnetic  permeability,  and  e  is  the  electric  permittivity. 

The  starting  point  for  the  completion  of  the  generalized  flux  function  in  Phase  II  is  the  assembled 
functions,  which  are  defined  below: 


dr 


dF  da  aj?'  , 


dF,  dG^  dHl' 
d^  dri  dC  _ 


+■7:21^2  +7731^3  +  ■nsSe  =  S, 


DNS 

LES 

CRF 

MHD 

CEM 

CAA 

ESI 

^l,i 

1 

0 

0 

0 

0 

0 

0 

^2,i 

0 

1 

0 

0 

0 

0 

0 

^3,i 

0 

0 

1 

0 

0 

0 

0 

^4,i 

0 

0 

0 

1 

0 

0 

0 

^5,i 

0 

0 

0 

1 

0 

0 

0 

^6,i 

0 

0 

0 

0 

1 

0 

0 

^7,i 

0 

0 

0 

0 

0 

1 

0 

^8,i 

0 

0 

0 

0 

0 

0 

1 

In  this  table,  the  subscript  “j”  on  n  has  been  introduced  to  identify  the  columns.  The  abbreviations 
in  the  table  are  as  follows: 


Symbol 

Description 

Remarks 

DNS 

Direct  numerical  simulation 

Standard  Navier-Stokes  Eqns 

LES 

Large  eddy  simulation 

Standard  Navier-Stokes  Eqns 

CRF 

Chemically-reacting  flow 

Continuum 

MOD 

Magnetogasdynamics 

Continuum 

CEM 

Computational  electromagnetics 
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CAA 

Computational  aeroacoustics 

FSI 

Fluid-structure  interaction 

Structure  equations  will  be 
included  in  the  Phase  II 
project 

The  complete  flux  functions  for  each  of  the  applications  are  subsets  of  the  following  generalized 
functions  derived  during  Phase  I: 


(tTi  +  7r2  +  TTg  +  7r4  +  TTe  +  7rT)(0  + 

(tTi  +  7r2  +  TTg  +  7r4  +  -TTs  +  7Vt)pU  + 

(tTi  +  7r2  +  -TTg  +  7r4  +  TTe  +  7rT)/?U  +  TTeB^ 

(tTi  +  7r2  +  -TTg  +  7r4  +  TTe  +  ^T)pW  +  TTg 

(tti  +  7r2  +  TTg  +  tit'jp  +  {^4  +  ^e)pBb  + 

-TVspYi  +  TTeDg  + 

(^4  + 


The  inviscid  fluxes  are 


/  (tti +  7r2 +  7rg +7r4  Ftte +  7r7)*or7  +  7r6(^^D^/e-4Z?y/e)  \ 

(tTi  +  7r2  +  TTg  +  ^j)  \puU  +  +  (^4  +  TTg) 

+7r6  ( ^  /  e  +  4  D  « /e) 

(tTi  +  7r2  +  -TTg  +  ‘nY)\pvU  +  +  (7r4  —  RjjB^^] 

(?ri  +  7r2  +7rg  +  T^Y)[pfwU  +  +  ('7r4  +  7rE)[4£&  - 

(tTi  +  ■7r2  +  -TTg  +  -TTy)  (/?£  +p)r7  +  (-^4  +  TT^)  ^EjjU  —  ■  B)^^] 

7rg[pyit7] 

(7^4  +  +  4^^)] 

\  (7r4  +  7rE)[(;i^,  +  u^j5,  +  / 


and 


/ 


V 


(tti  +  tt 2  +  TTg  +  -7:4  +  -tte  +  Trr)^^  +  -TJ^DylY} 

(tTi  +  7:2  +  TTg  +  TTt)  [^uV  +  7]  Jj]  +  (7:4  +  TTe) 

+7r6(-7],Z),/e  +  7],Z)Je) 

(tTi  +  772  +  TTg  +  ^t)  +  7]^;??]  +  (774  +  7^B)[%Eh  “ 

+776(77,Dje-7j^D,/e) 

(771  +  772  +  77g  +  77T)[pii^  V  +  7]^p]  +  (774  +  77E)[7]^£i  -  RjhBy 

(771  +  772  +  77g  +  77t) (p-E  +  p)  V  +  (774  +  77e) [2£i V  -  (U  ■  B) 
+776 

77g[/7l-V]  +  776  (-7J^B^//^^ 

(774  +  77E)[u7Jy  +  -  u{7]^By  +7]^B^)] 

(774  +  77s)[(u7^^  +  W7].^)By  -  V  (7]^B^  +  7]^B^)] 

(^4  +  77E)[(7i77,  +  ;;77jB,  -^i;(77,B, +7^^^^)] 


\ 


y 
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/  (tti +  7r2 +7r3 +7r4  +  7rs +  7r7)i0iy +  7r6(Cj^Z:)^/e- \ 

(tti  +7r2  +7r3  +  ^t)[puW  +  Cjj]  +  + 

+^6(-C^s/^+ A) 

(tTi  +  7r2  +  TTa  +  TTy)  [puiy  +  CyP]  +  ('^4  +  ■TTs)  [Cy^b  - 

+ne{C^Dy/e -CyD^^.M 

(tti  +  TTa  +  tts  +  +  Cj?]  +  (?r4  +  [C^i  - 

+7r6(-Cj;^^/Ai^  + 

(tti  +  7r2  +  tts  +  7r7)(p^  +^j)Ty  +  (7:4  +  Tr^)[2Ej,W  -  i?i(U  ■  B)^] 

n3\pyiW]  +  ne(-(:^By/^^  +  CyB^/^^) 

(7r4  +  7r£)[uC5,  +wCjB^  -  uiC^B^  +  C^B^)] 
iTT4  +  7r^)[{uC^  +wCJBy  -v{<;^B^  +  CA)] 

\  {n4  +  ng)[luC^+vCy)B^ -w(C^B^  +  CyBy)]  j 


The  dissipative  fluxes  are  as  follows: 


[(’■‘1  +  ''r2  +  TTs  +  7:4  +  TTe  )Til  +  (772  +  TTg  )c7fi®®] 

[(■'Tl  +  ■'''2  +  TTs  +  774  +  TTe  )Ti2  +  (772  +  77a 
m  [(771  +  772  +  773  +  774  +  775)  Tia  +  (772  +  773) 

[(77l+  772  +  773  +  774  +  77E)rij  +  (772  +  773)tTfj-®®] 

+  (Tl)pij7fH.^^.  [(’Tl  +  772  +  773  +  774  +  77e)®  +  (772  +  773)f2f ^®] 

+K  +  77e)^L.^(|^)-4..(*^  +  4^+£^)] 


K  +  77e)^ 
(’^4  +  775):^^ 


y  ^ 

e  9  B. 
3x 


9iS  /. 


+  l 


3  B. 


3  gm 

■=^3&  ju„ 


’s;  3si  ij,^ 

Se 


-4 


-hVoy,  [(^1  +  ?r2  +  ^Ts  +  7r4  +  7rs)Tii  +  {^2  + 

[(■^1  +  ^^2  +  '^S  +  ■7r4  +  '^s)T'i2  +  ('^2  +  -TTS ) ] 

[('^1  +  7r2  +  TTg  +  7:4  +  7rE)Ti3  +  (7:3  + 

1^^==,  [(^1  +  7r2  +  TTg  +  7:4  +  ns)Ti^  +  (7:3  + 

[('^1  +  7r2  +  TTa  +  7:4  +  TTs)  $  +  (7:2  +  TTa  )(3f 


^(■7-l)Pr  M=  Ri 

+  (7:4  +7r£)- 


TT  I  /aJ-f  3k;(  \2  )  '>Xi\  5s  ^  3y  '  3s  J\ 

Pij  +  %f/^)] 


(^4  +  7rs)^ 
(7r4  +  7rE)^ 
('^4+TrE);^ 


V 

%3y/^, 

V 

V 


.  J_  ^  _^:^7E _ ^ 


+  77 
^  Vl/  3y 


“"^45!- 


■-^4fe 

■  -77 

/y  Se  i+^ 


and 
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K 


( 


1 

j 


V 


0 

[K  +  ^2  +  TTg  +7r4  +7rs)Tii  + 

[(^1  +  TTs  +  TTg  +  7r4  +  Tr^)ri2  +  (^Ts  +  ^3)cf2^^] 
r;C^  [i^l  +  TTs  +  TTg  +  7r4  +  7rs)Ti3  +  (^2  + 

^  (C^  U,  [(tTi  +  TTs  +  TTg  +  7r4  +  7rE)Tij  +  (tTs  + 

+  (7-1)  rSs!,  [(■^l  +-7r2  +'^3  +  '^4  +  -^s)^  +  ('^2  +'^3)Qi 

^  SW 

^  I  ^  ^ 

3x  ^ 

^  S  _|_  /*  _^  ^  S 


(^4+^e):^ 
(7r4  +  7rE)^ 


\ 


/ 


The  proposed  software  suite  involves  the  development  of  a  systematic  way  of  1)  computing 
derivatives,  2)  carrying  out  TDMA,  pentadiagonal  solution  (PDMA),  or  Jacobi/Gauss-Seidel 
iterations,  and  3)  updating  the  residuals.  The  schematic  shown  below  represents  the  initial 
concept  for  the  multidisciplinary  software,  while  Appendix  A  shows  a  detailed  design.  The 
following  remarks  are  in  order.  Note  that  within  the  “engine,”  each  application  (CEM,  MHD, 

ESI,  CAA,  etc.)  could  have,  as  an  alternative  to  the  “tt”  formulations  above,  a  separate  routine 
that  accesses  the  same  derivative  calculation  function.  The  Phase  I  project  focused  on  the  “tt” 
formulation.  An  alternative,  optimized  approach  will  be  explored  in  Phase  II  in  order  to  obtain  the 
most  efficient  code.  This  is  intended  for  the  applications  that  have  less  number  of  variables  than 
the  8  in  the  base  code.  The  graphics  user  interface  (GUI)  provides  the  appropriate  software 
branches,  as  is  the  case  in  the  current  version  of  the  commercial  code  INSTED.  We  will  also  have 
core  and  application-specific  results  and  visualization  (post-processing). 

Programming  in  Phase  I  was  based  on  FORTRAN  90,  C,  and  C++.  The  preliminary  compute 
engine  used  FORTRAN,  whereas  C  and  C++  were  used  mostly  for  the  interface  codes. 

IV:  Graphical  User  Interface  (GUI) 

Although  not  promised  for  Phase  I,  we  implemented  a  preliminary  graphics  user  interface  (GUI), 
to  demonstrate  the  feasibility  of  a  user-friendly  interface  for  the  proposed  multidisciplinary 
software.  This  capability  is  an  absolute  necessity  for  commercialization  of  the  final  product.  The 
GUI  contains  most  aspects  of  the  proposed  software,  except  the  engine  (solver).  That  is,  the  left- 
hand  side  of  the  schematic  shown  below.  Sample  screens  applicable  to  one  of  the  applications  are 
shown  in  Appendix  B. 
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Figure  1  Preliminary  Multidisciplinary  Software  Concept 
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V.  Design  of  the  multidisciplinary  Software  Suite 

A  detailed  design  of  the  multidisciplinary  software  is  shown  in  Appendix  A. 

VI.  Demonstration  of  the  multidisciplinary  Software  Suite 

The  feasibility  of  the  proposed  multidisciplinary  software  suite  was  established  in  Phase  I,  by 
ascertaining  accuracy  and  computational  efficiency  when  the  proposed  design  in  Appendix  A  is 
implemented.  The  base  code  for  the  development  of  the  multidisciplinary  code  is  the  Air  Force 
code  FDL3D/MGD,  which  has  a  default  number  of  variables  of  8.  We  have  modified  this  code  to 
execute  in  a  multidisciplinary  environment  and  have  used  two  problems  for  evaluation: 

a)  MGD  (Alfven  wave  with  ohmic  damping).  Here,  we  compare  the  computational  efficiency 
and  accuracy  of  the  original  FDL3D/MGD  code  with  those  of  the  preliminary 
multidisciplinary  code  derived  from  it.  The  CPU  times  are  shown  in  Table  I  and  the  accuracy 
below  in  figures  2  through  4.  From  the  table  we  see  that  the  efficiency  of  the 
multidisciplinary  code  is  comparable  to  that  of  the  original  MGD  code  for  both  the  Runge- 
Kutta  and  Beam-Warming  time  integration  methods.  Note  that  the  CPU  times  shown  are  for 
100  time  steps  and  the  Beam- Warming  calculations  use  3  sub-iterations  per  time  step.  This 
exercise  shows  that  the  multidisciplinary  software  suite  is  quite  feasible,  in  relation  to  the 
base  code. 

b)  Acoustic  scattering  by  a  concept  X24C  vehicle.  This  problem  is  interesting  for  a  couple  of 
reasons:  it  represents  a  complex  geometry  and  involves  only  5  variables  (as  opposed  to  the  8 
in  the  base  code  for  the  multidisciplinary  software).  Therefore,  the  associated  overhead  is 
expected  to  challenge  the  multidisciplinary  code  relative  to  a  dedicated  code  for  acoustic 
scattering.  (See  Ladeinde  et  ah,  2001.)  The  initial  results  confirm  our  suspicion  (not  shown). 
However,  a  preliminary  optimization  technique  was  found  to  show  great  promise  for  the 
multidisciplinary  software  suite  when  the  number  of  variables  is  much  fewer  than  that  in  the 
base  code.  The  latter  will  be  fully  developed  in  Phase  II. 
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Lilts  from  Multi-disciplinary  code  using  Runge-Kutta  scheme 
Mesh  System:  101x3x3) 


its  from  Multi-Disciplinary  code  with  B 
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APPENDIX  A 


Design  of  the  Multi-disciplinary  Software  Suite 


•  New  Control-Variable  Table 


DNS/NS 

LES/NS 

Combustion 

CAA 

DNS/MHD 

LES/MHD 

CEM 

ni 

1 

m 

1 

n3 

1 

U4 

1 

U5 

1 

U6 

1 

m 

1 

•  Flow-Chart 
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•  Modification  on  Highlighted  Blocks 


Modifications: 

1.  In  “Gconinion.f90”  add  code  segment: 

I 

I 

I 

I 

I 

INTEGER: :  MULTIDSCP 

INTEGER,  DIMENSION(7) ::  MULTI  TYPE 

I 

I 

I 

I 

I 

2.  In  “I 


MULTI_TYPE=0 
READ  *,  MULTI  DSCP 
MULTI_DSCP=MULTI_DSCP/1 0 
MULTI  TYPE  (MULTI_DSCP)=1 
WRITE  *,  “WARNING  SIGNS” 
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•  Modification  on  Highlighted  Blocks 


20 


DOD/SBIR  Phase  I  Final  Report 

Advanced  Software  Suite  for  Multidisciplinary  Computation  (Thaerocomp  Technical  Corp.) 


March  3,  2003 


•  Modifications  on  Transformation  between  Primitives  and  Conservatives 


PV  to  CV 


CALL  CONVERT 
Or  CALL  DGCONVERT 

Or  CALL  UPDTQO 


CV  to  PV 


> 


CALL  UDGCONVERT 


COMMON  CONTROL  LOGICS: 


VARIABLE_ARRAY=0  !  INITIALIZATION 

! 

!  CASE-SELECT  STRUCTURE 
! 

MULITI  ONE:  SELECT  CASE  (MULTI  DSCP) 

CASE  (1:6) 

VARIABLE_ARRAY=BASIC  TERMS  in  NS 
CASE  (7) 

VARIABLE_ARRAY=CEM  TERMS 
END  SELECT  MULTI  ONE 

MULITI  TWO:  SELECT  CASE  (MULTI  DSCP) 

CASE  (3) 

VARIABLE_ARRAY=REACTING  SPECIES 
CASE  (4:5) 

VARIABLE_ARRAY=B  FIELDS 
END  SELECT  MULTI  TWO 


Notes: 

1.  Replace  the  VARIABLE  ARRAY  by  different  WORK  ARRAYS  in  different  subroutines, 
namely, 

a)  TEMP  0  in  CONVERT.F90; 

b)  RHS  0  in  DGCONVERT.F90 

c)  QOLD  0  in  UPDTQO.F90 
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•  Modifications  on  Beam-Warming  Solver 
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CALL  SWEEPJl 

!  (CALL  RHSXZ) 

CALL  SWEEPI 

!  BY  IVMON 

CALL  SWEEPK 

!  (CALL  RHSY) 

CALL  SWEEP J2 
(CALL  SDAMP/X/Y/Z/  in 

SWEEPK) 

COMPACT: 

CALL  VCMPYRHS 
CALL  VCMPXZRHS 

2"'*  ORDER: 

CALL  V2NDXZRHS 
CALL  V2NDYRHS 
CALL  VCROSX 

CALL  VCROSY 

CALL  VCROSZ 
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4. 


CALL  CMATA 
CALL  CMATB 
CALL  CMATC 


5. 


CALL  CMATVl 

dFy  dGy  dHy 

CALL  CMATV2 

< - > 

CALL  CMATV3 

dQ,  ’  dQ^  ’  dQ^ 

6.  The  same  CASE-SELECT  logics  as  used  above  are  applied  here.  The  work  arrays  specifically  are: 

a)  WRK  (*,*,*,*1),  TEMP  (*,*,*)  and  RHS  (*,*,*)  in  SWEEP/1/J1/J2/K/.F90 

b)  WRK  (*,*,*,1)  in  /X1/ET/ZT/COMP.F90 

c)  WRK  (*,*,*,1)  and  RHS  (*,*,*)  in  VCMP/XZ/Y/RHS.F90  and  VCROS/X/Y/Z/.F90 

d)  ALCL  (*,*,*)  in  CMAT/A/B/C/.F90  and  CMATW1W2/V3/.F90 

e)  WRK  (*,*,*,*)  in  SDAMP/X/Y/Z/.F90 

7.  Note:  Roe,  Upwind,  2’'‘*-Order  differencing  schemes  are  currently  not  implemented  for  the  multi¬ 
disciplinary  computation. 

8.  ADD  SUBROUTINE  MULT1_SOURCE.F90  for  Beam-Warming  procedure.  The  algorithm  is: 


LHS  ofBW  (0) 

At 

=RHS  of  BW  (0)  + - 

l-hO 


Where  BW  (0)  represents  the  original  Beam- Warming  procedure  and  the 
expression  of  the  source  term  S  is  described  in  my  notebook. 
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•  Modifications  on  Runge-Kutta  Solver 


Detail  Descriptions  of  the  implementation: 

1.  . . : 


Pre-Processing 


Post-Processing 


2. 


CALL  RUNGE 
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3.  The  same  CASE-SELECT  logics  as  used  above  are  applied  here.  The  work  arrays,  specifically,  are 

a)  TEMP  (*,*,*)  and  QOLD  (*,*,*,*)  in  RUNGE.F90 

b)  WRK  (*,*,*,1)  in  /XI/ET/ZT/COMP.F90 

c)  WRK  (*,*,*,1)  and  RHS  (*,*,*)  in  VCMP/XZ/Y/RHS.F90  and  VCROS/X/Y/Z/.F90 


4.  ADD  SUBROUTINE  MULTI_SOURCE.F90  for  Runge-Kutta  procedure.  The  algorithm  is: 


LHS  ofRK(O) 

=RHS  of  RK  (0)+  At- J-S 

Where  RK  (0)  represents  the  original  Runge-Kutta  procedure  and  the 
expression  of  the  source  term  S  is  described  in  my  notebook. 
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Appendix  B:  Sample  GUI  Screens 
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Table  I:  Comparison  of  CPU  time  (seconds)  for  the  calculation  of  Alfven  wave  with  ohmic  damping.  The 
CPU  times  are  reported  for  100  time  steps.  The  Beam-Warming  scheme  uses  3  sub-iterations  per  time  step. 


Code  Type 

Grid  Size 

CPU  time  (seconds) 

Runge-Kutta 

Beam- Warming 

Base  code 

101x11x11 

3035.9 

4742.1 

Multidisciplinary  code 

101x11x11 

2990.9 

4787.0 

Multidisciplinary  code 

101x3x3 

3632.0 

95.09 
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Appendix  C:  List  of  FDL3D  Files  Modified 

The  following  is  a  list  of  FDL3d/MGD  project  files  that  were  modified  in  the  preliminary 
conversion  of  the  code  to  a  multidisciplinary  suite.  Note  that  some  files  contain  multiple 
routines/procedures. 

bndry.f90,  cmata.f90,  cmatb.f90,  cmatc.f90,  cmatvl.f90,  cmatv2.f90,  cmatv3.f90,  convrt.f90, 
datin.f90,  deriv.f90,  dgconvrt.f90,  etcomp.f90,  fdl3d_ice.f90,  fdtry.f90,  gcomonl.f90, 
met2nd.f90,  multi_src.f90,  nirhs.f90,  output.f90,  pirhs.f90,  rhsrunge.f90,  runge.f90, 
sdampx.f90,  sdampy.f90,  sdampz.f90,  sdampzl.f90,  sdiagx.f90,  sdiagy.f90,  sdiagz.f90, 
sweepi.f90,  sweepjl.f90,  sweepj2.f90,  sweepk.f90,  tmstep.f90,  txiirhs.f90,  tzetrhs.f90, 
udgconvrt.f90,  updtqo.f90,  vcmpxzrhs.f90,  vcmpyrhs.f90,  writetec.f90,  wrtfield.f90, 
wrthist.f90,  xicomp.f90,  ztcomp.f90 
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APPENDIX  D:  Grid  for  Second  Test  Problem  for  the  Multidisciplinary  Software 


Fig.  5  Coarse  grid  model  of  X24C  Reentry  Vehiele  showing  the  surfaee  mesh. 
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